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Abstract. We investigate dynamics of a flat FRW cosmological model with a barotropic matter and 
a non- minimally coupled scalar field (both canonical and phantom) . In our approach we do not assume 
any specific form of a potential function for the scalar field and we are looking for generic scenarios of 
evolution. We show that dynamics of universe can be reduced to a 3-dimensional dynamical system. 
We have found the set of fixed points and established their character. These critical points represent 
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(rapid- roll and slow- roll) , (c) a radiation domination, (d) a matter domination and (e) a quintessence 
era. We have shown that the inflation, the radiation and matter domination epochs are transient ones 
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effect of a non-minimal coupling constant. We show the existence of a twister type solution wandering 
between all these critical points. 
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1 Introduction 

In modern cosmology a scalar field <fi plays a very special role. The discovery of accelerated cosmic 
expansion [1, 2] gave a motivation to study dynamics of dark energy models (see [3] for review). In 
the context of the quintessence idea [4, 5] the simplest dynamical models involving the scalar field <f> 
with a potential function V(<f>) are used to model a time dependent equation of state parameter w^. 
While the simplest candidate for dark energy seems to be a positive cosmological constant, the ACDM 
model is favoured by observational data [6-8], such a explanation of cosmic acceleration suffers from 
the fine tuning problem [9] and the coincidence problem [10]. In order to alleviate those problems 
many alternatives have been proposed like phantom dark energy [11, 12] or extended quintessence 
[13-16]. 

If we are going to generalise the scalar field cosmology minimally coupled with gravity, then 
inclusion of the non-minimal coupling term of type —^R<p 2 [17-19] seems to be natural and the 
simplest generalisation of the Lagrangian for scalar field dynamics in the background of cosmological 
models with maximal symmetry of space-like slices. Of course the value of this additional parameter 
should be estimated from observational data [20-22] or given from some theoretical arguments [14]. 
The nonzero £ arises from quantum corrections [23] and it is required by the renormalization [18]. 
While the simplest minimally coupled scalar field with a quadratic potential function has strong 
motivations in observations [24, 25] its generalisations with a non-minimal coupling term was studied 
[26] in the context of origin of the canonical inflaton itself. 

In this paper we investigate the dynamical evolution of scalar field cosmological models with 
a non-vanishing coupling constant between a scalar field and gravity. The role of the non-minimal 
coupling in evolution of the universe in the context of inflation and quintessence was studied previously 
by many authors [27-48] and in connection with the development of the Standard Model with a non- 
minimally coupled Higgs field [49-52]. The dynamical systems methods are used in investigating of 
evolutional paths of cosmological models which dynamics is parameterised by the energetic variables 
very useful in this context [53]. 

We are adopting dynamical systems methods in studying the evolution of the cosmological model 
and its dynamics can be visualised in the phase space which is a geometric framework for its explo- 
ration. Moreover one can investigate all evolutional paths of the system under consideration for all 
admissible initial conditions. Therefore one should ask whether different models with a desired prop- 
erty are typical (generic) in the class of all possible models. In our opinion physically interesting 
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cosmological models should be generic in the sense that they do not depend on the special choice of 
initial conditions which should be determined from quantum models. 

To keep generality of our considerations we do not assume any specific form of the potential 
function of a canonical or phantom scalar field. It will be demonstrated that the parameter of 
non-minimal coupling £ plays the crucial role during the cosmological evolution. We will show the 
emergence of a new phase space structure organised through critical points and trajectories. For 
completeness we also include the model with the barotropic matter with the constant equation of 
state parameter w m . 

We demonstrate that, in principle, the phase space structure at a finite domain is determined 
by five critical points corresponding to important events during the cosmic evolution, namely, the 
singularity (of a finite scale factor type), inflation, radiation and matter dominated epochs and finally 
the accelerated expansion era. In our previous paper we introduced notion of the twister solutions 
the solutions linking the subsequent cosmological epochs [54]. In the present paper we generalise 
this notion without assuming any form of the potential function of the scalar field. The evolutional 
scenarios investigated in this paper are obvious only if the non-minimal coupling is different from 
minimal (£ = 0) and conformal (£ = 1/6) coupling value. In this sense we study a unique type of 
evolution. 



2 The model 

In the model under consideration we assume the spatially flat Friedmann-Robertson- Walker (FEW) 
universe filled with the non-minimally coupled scalar field and barotropic fluid with the equation of 
the state coefficient w m . The action assumes following form 

= 11 d A x^l^R-e(g^d^d^ + ZR^-2U(<j ) )]+S m , (2.1) 



2 

where k 2 = 8irG, e = +1, —1 corresponds to canonical and phantom scalar field, respectively, the 
metric signature is (—,+,+,+), R = 6^ + is the Ricci scalar, a is the scale factor and a dot 

denotes differentiation with respect to the cosmological time and U((j>) is the scalar field potential 
function. S m is the action for the barotropic matter part. 

The dynamical equation for the scalar field we can obtain from the variation 6S/64> = 

4> + 2>H4> + £R4> + eU'(<f>) = 0, (2.2) 

and energy conservation condition from the variation SS/Sg^ = 

£ = £;U 2 + e3£H 2 (j) 2 + sZ£,H{4> 2 ) + U(<f>) + p m - -\rH 2 . (2.3) 

2 K z 

Then conservation conditions read 

-%H 2 = p^ + p mi (2.4) 



H = - — 
2 



{P(j> + P<t>) + Pm(l + W m ) (2.5) 



where the energy density and the pressure of the scalar field are 

P4> = s\ft + U{4>) + e3(H 2 cj) 2 + e^H^ 2 ), (2.6) 
P* = - - U M + ^H(^ 2 )--e2ai - H)H<f - e3£(l - 8£)#V + 2&U' (</>). (2.7) 

Note that, when the non-minimal coupling is present, the energy density and the pressure p^, 
of the scalar field can be defined in several possible inequivalent ways. This corresponds to different 
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ways of writing the field equations. In the case adopted here the energy momentum tensor of the 
scalar field is covariantly conserved, which may not be true for other choices of and [55, 56]. For 
example, the redefinition of gravitational constant = k~ 2 — e^4> 2 makes it time dependent. The 
effective gravitational constant can diverge for a critical value of the scalar field <fi c = ±(eK 2 ^) -1 / 2 . 
Though the FRW model remains regular at this point, the model is unstable with respect to arbitrary 
small anisotropic and inhomogencous perturbations which become infinite there. This results in the 
formation of a strong curvature singularity prohibiting a transition to the region k 2 s < [57]. 
In what follows we introduce the energy phase space variables 



which arc suggested by the conservation condition 



(2- 



flj, + O m — 1 



or in terms of the newly introduced variables 

ft = y 2 + e [(1 - 6£)x 2 + 6£(x + z) 2 
The acceleration equation can be rewritten to the form 



= 1-0; 



H = —[PeS +PcS) = (1 + ffi eff ) 



(2.9) 



(2.10) 



(2.11) 



where the effective equation of the state parameter reads 
1 



WcS 



1 - e6f (1 - 6£)z 2 



1 + e(l - 60(1 - w m )x 2 + e2f(l - 3w m )(x + zf + 
+ (1 + w m )(l - y 2 ) - e2£(l - 6£)z 2 - 2£y 2 \z 



(2.12) 



where A - ^ l duj£ 

K U{4>) d0 ' 



The dynamical system describing the investigated models is in the following form [54, 58] 



r6£(l-6£)z 2 +_(a. + 6fxt) 



3 - 2CA2/ 2 z 



+e(l - 6^(1 - w m )x 2 + e2£(l - 3w m ) (x + zf + (1 + w m ){l - y 2 ) 



2/(2- -Xx 



1 - e6£(l - 6£)z S 



2£Xy 2 z 



A' 



+e(l - 60(1 - w m )x 2 + E 2£(l - 3w m ) (x + zf + (1 + w m )(l - y 2 ) 

1 -£60l-60z 2 
A 2 {T -l)x\l- e^{l-Q£)z 2 



where a prime denotes differentiation with respect to time r defined as 



_d 

d7 



1 - e6f (1 - 6£)z' 



dln< 



(2.13a) 

(2.13b) 

(2.13c) 
(2.13d) 

(2.14) 



where the expression in brackets is assumed as a positive quantity to assure that during the evolution 
with r > the scale factor a is growing, i.e. the universe expands, and 

U"{cp)U{<t>) 

U'{4>) 2 ' 
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parameters 



z(X) 



potential function U(4>) 



t/ exp (-f<£ 2 
C/ o exp(^exp(/30) 



a^O, ,3 = 0, 7 = 
a = 0, 0^0, 7 = 

a^Q, /3^0, 7 = 

a^O, /3 = 0, 77^0 
a = 0,^0,7^0 

a^O,/?^ 0,7^0 



— + const. 



In A 

P 



const. 



const. 



]^+m+ const 
'(^) 



arctan ( 

3 + 2^A 



const, 
const. 



2 arctan 



V-)j2 + 4 Q 



x /- f 3^+4a 7 



const. 



- const. <j 



>W)) 



1 

Uo (7^ — const.)"' 
Uq exp (i (a0 + conts. exp ( 

Uq (cos (y/orf((f> — const.))) 7 
(exp (const./?) + 7 exp (/?</>)) T 

(l^) ( cos (i\/-/3 2 +4a 7 (0 - const. 



7) 
?7o 



C/o exp 



Table 1. Different examples of potential functions for various configurations of parameters values of the 
assumed form of the F(A) function T(A) = 1 — -A- (a + f3\ + 7A 2 ). 



To investigate the dynamics of the universe described by the dynamical system (2.13) we need to 
define an unknown function T, i.e. we need to define the potential function U(4>). In the special 
cases of the system with the cosmological constant or exponential potential, U = Uq = const, or 
U = Uq exp (— A</>), the dynamical system (2.13) can be reduced to the 3-dimensional one due to the 
relation A = and T = in the former case, and A = const, and T = 1 in the latter case. Then 
dynamical system consists of three equations (2.13a, 2.13b, 2.13c). 

There is another possibility of reduction of the system (2.13) from a 4-dimensional dynamical 
system to a 3-dimensional one. If we assume that z = z(X) and T = T(A), then using (2.13c) and 
(2.13d) we can find the function z(A) from the differential equation 



dz(A) 



z'(X) 



1 



dA v ' A 2 (r(A)-l 

which can be integrated for some given function T(A) 

z(X) = - 



(2.15) 



dA 

A 2 (r(A)-i)' 

Then the dynamical system describing the investigated models is in the following form [54, 58] 



(2.16) 



> - e\\y 2 ) [l - e6f (1 - 6£)z(A) 2 ] + 1 ( x + 6£z(A)) 



+e(l - 60(1 - w m )x 2 + e2£(l - 3w m ) (x + z{\)f + (1 + m m )(l - y 2 ) 



2 - X -\x 



1 - e6£(l - §t)z(X)< 



2 V 



- - 2£\y 2 z(\) 



A' 



+e(l - 60(1 - w m )x 2 + e2£(l - 3w m ) {x + z(A)) 2 + (1 + m m )(l - y 2 ) 

-A 2 (r(A)-l)x[l-e6£(l-6£MA) 2 

where a prime denotes now differentiation with respect to time r defined as 

d 



_d_ 

d7 



1 - e6£(l - 6£)z(A)' 



din a 



(2.17a) 

(2.17b) 
(2.17c) 

(2.18) 



and we assume that the term in bracket is positive in the phase space during the evolution. 
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Now we are able to express the acceleration equation (2.11) in terms of the energy phase space 
variables and time r 

= -3 [l - e6f (1 - 6£XA(t)) 2 ] (1 + w cS ) (2.19) 
which together with (2.12) results in 



In (t^t) 2 = -3 /; {l + w m + e(l - 60(1 - w m )x(r) 2 + e2£(l - 3«7 m )(x(r) + z(A(t)))' 
-2,(t) 2 (2£A(t)z(A(t)) + 1 + w m ) - e8£(l - 60^(A(r)) 2 }dT, 



(2.20) 



where i? ml denotes the initial value of Hubblc's function at time r = 0. In what follows we will be 
using this expression together with the linearised solutions in the vicinity of every critical point to 
investigate the behaviour of Hubblc's function with respect to the scale factor, as well as Hubble's 
radius defined as 

Rh = h- 

For example if the function T(A) is assumed in the following form 

r(A) = l-^(a + /3A + 7 A 2 ), 

then in Table 1 we have gathered forms of the functions z(X) and corresponding potential functions 
U{(j>) for various configurations of values of parameters a, /3 and 7. As we see there are various 
potential functions which are the most common used in the literature of the subject. Of course this 
simple ansatz for the function T(A) does not manage all possible potential functions. Let us consider 
the following function 

r W -5 ™ , 

i 4(2 ± Vi + v 2 \ 2 ) 

as one can check from (2.16) we receive 



2 ± V4 + cx 2 A 2 

2(A) = h const. 

A 

and this example corresponds to the Higgs potential 

U(4>) = TjM<f>- const.) 2 -cr 2 

We need to stress that the discussion presented below is not restricted to the specific potential function 
but is generic in the sense that it is valid for any function T(A) for which the integral defined in (2.16) 
exists. 



3 Dynamics of universe with a potential 

In our investigations of dynamics of the system given by eqs. (2.17) we will restrict ourselves to the 
finite region of the phase space, i.e. we will be interested only in the critical points in the finite 
domain of the phase space. The full investigations of the dynamics requires examination of critical 
points at infinity, i.e. compactification of the phase space with the Poincare sphere. The procedure 
of transforming dynamical variables into the projective variables associated with the compactification 
requires that the right hand sides of the dynamical system should be polynomial. In our case this 
is not always true because of the form of function z(\) (see table 1). In what follows we present 
detailed discussion of character of critical points of the system (2.17) corresponding to different stages 
of cosmological evolution (table 2). 
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X* 


y* 


A 


Wcff 


1. 


^i = -6£2(A X ) 




\* "\ \2 1 


±00 


ZcL. 

2b. 


x 2a — ~ bt i z { A 2a) 


(-,.* ^2 _ 45 
\i)2a) — 2fA* z(A* ) + (l+wi m ) 
f * \2 _ ' 2£(l-3u> m ) 


\* . .,( \\2 _ 1 
A 2a ■ — e6f(l-60 


in AC 

w m -2(, 


(l-6«)(2{A5^(A| ! ,)+(l+«m)) 


^26- 2 (^) — e 64(l-64) 


1-64 


3a. 
3b. 


4a = 9(x) = 1 
^36 = 


l/3« = o 

Z/3& = 


^3a ■ Z W — £64(1-64) 
^3b : Z W = £64 


1 

3 
1 
3 


4. 


x* 4 = 


wS = o 


A4: 2 (A) = 


W m 


5. 


x* 5 = 


(^) 2 = l- £ 6^(A^) 2 


A^: \z(\y+4z(\)-^=0 


-1 



1 g{x) = e(l - 44 - w m )x 2 + e44(1 - 3w m )z(\* a )x + T ||j(l - 3w m ) 



Table 2. Critical points of the system under consideration. 
3.1 Finite scale factor initial singularity 

Our discussion of the dynamics of the model under consideration we begin with the critical point 
located at 

xt = -6^z(\t), yt=0, X{: z(Xf = -^-^ (3.1) 

where the last expression means that the coordinate A^ of the critical point is the solution to the 
equation z(X) 2 = jg^r^g^y- This critical point represents a singularity because the value of w e g given 
by (2.12) calculated at this point is 

W e f[ = ±00. 

We need to stress that this critical point exists only if e£(l — 6£) > 0, i.e. for the canonical scalar 
field (e = +1) for < £ < 1/6, and for the phantom scalar field (e = —1) for £ < or £ > 1/6. 

In cosmological investigations one encounters usually various types of singularities such as: ini- 
tial finite scale factor singularity [59, 60], future finite scale factor singularities: the sudden future 
singularities [61], the Big Brake singularity [62], and the Big Boost singularity [63]. 

Linearised solutions in the vicinity of this critical point are 

an(r) = x\ + ((zf - ^) +2(1 - 3£)z'(A 1 )(Af - A?)) cxp (l lT ) 

-2(1 - 3f)/(Ai)(Af - XI) cxp for), 
yx(r) = j/rexp(/ 2 r), 
Ai(r) = A* + (Ar-A* 1 )cxp(/ 3 T). 

where 

h = 6f, h = 6£, h = 12£ 

are the eigenvalues of the linearization matrix calculated at this critical point, x™ 1 , y™ 1 and A 1 / 11 are 
initial conditions. For positive values of the coupling constant £ > this critical point represents 
an unstable node type critical point. For £ < which is possible only for the phantom scalar field 
(e = —1) the critical point is of a stable node type. 
Using the time reparameterization (2.18) 

dlna = (l - e6£(l - 6£)z(A(t)) 2 ) dr (3.3) 

and expansion into the Taylor series around the critical point coordinate A* 

z(A) = z(A*) + z'(A*)(A-A*) 

up to linear terms 

z(X) 2 = z{X*) 2 + 2z{X*)z'(X*){X - A*) 



(3.2a) 
(3.2b) 
(3.2c) 
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and then together with (3.2c) we have 

z(X(r)f = z(X* 1 ) 2 + 2z(X* 1 )z'(X* 1 )(X'f - A^expfor). 
Inserting this expansion into the time reparameterization we receive 

dlna = -el2f(l - 6£)z(\l )z , (X* 1 )(X™ - A*) exp (/ 3 r)dr 
which can be directly integrated for Z 3 = 12£ > 

,o 

Alna = -el2e(l-602(A^)z'(A*)(A I 1 m -A*) / exp (Z 3 r)dT 

J —oo 

(we could take also Z 3 = 12£ < and the integration in this expression should be taken (0, 00) because 
for £ < this critical point represents a stable node). Where the result is 

/o ini \ 

Alna = ln( -M = — e(l — 6^)z(A 1 )z'(A 1 )(A i 1 ni — A*) 

Finally we receive 

< = a s exp { - e(l - BZ)z{K)z , (\l)(X? i - A*)} 

where the value in the exponent is finite, and a™ 1 is the value of the scale factor at r = and a s is the 
value of the scale factor at singularity. This equation gives us the scale factor growth from singularity 
to the some initial point where linear approximation is still valid. From simple considerations we have 
that 

z(Xl)z'(Xl)(\r-Xl)<0 

which show that the critical point under consideration represents the finite scale factor singularity. 
Moreover it is a past singularity for the canonical scalar field with < £ < 1/6 and the phantom 
scalar field with £ > 1/6 and future singularity for the phantom scalar field with £ < 0. 

Now, using linearised solutions (3.2), we can express (2.18) and (2.20) as parametric functions 
of time t 

f In (^) = s(l - 6£>(A!>'(A!)(Af » - XI) (I - exp (Z 3 r)), 

J In ( 7 ^) 2 = 3{ - 4£r - 4(1 - eOziXDl^f - *I) + 2(1 - S^TO^f - A|)](l - exp&r)) 
[ +4(1 - 60(1 - 3w m - 12£)z(At)z'(A 1 )(Af - A|) (l - exp (f 3 r))}. 

(3.4) 

The linearised solutions used to obtain these relations are valid up to the Lyapunov characteristic 
time which is equal to the inverse of the largest eigenvalue of the linearization matrix. In our case it 
is T cnc i = 7^ = Inserting this in to equations (3.4) we obtain maximal values of the scale factor 
and the Hubble's function respectively: 

'ln(^) = e (l-6£)z(A^TO(Ar-At)(l-e), 

< hi (fw) 2 = 3{ - i - £ |(1 - 6£MAT)[(zf - O + 2(1 - 3£K(A|)(Ar - A|)](l - e) (3.5) 
+4(1 - 6£)(1 - 3w m - 120z(X* 1 )z / (X* 1 ){Xf i - AJ)(1 - e)}. 

The plot representing the evolution of these quantities together with Hubble's horizon is presented in 
figure 1. As one can simply conclude 

£ > : lim H' 2 -> 00. 

T— >- — OO 

The general conclusion is that any phantom scalar field cosmological model with the negative 
coupling constant £ < and the potential function which can be represented by function z(X) possesses 
the finite scale factor future singularity with w c s = ±00. 
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Figure 1. Evolution of InH 2 (left panel) and Rh (right panel) as a function of a natural logarithm of the 
scale factor lna for a sample trajectory with e = +1, £ = |, z'{\\) = — = 100 in the vicinity of the critical 
point corresponding to the finite scale factor singularity. The solid black line represents the linearised solution 
(3.4) and the dotted line corresponds to the numerical solution of the system (2.17). 



3.2 Inflation with the non-minimal coupling and arbitrary potential 

Now we proceed to the very important phase of the evolution of the universe, namely the inflation. 

3.2.1 Fast-roll inflation 

The critical point located at 

4£ 



2a 



2^X* 2a z(X* a ) + (l + w m )' 



A 



2a 



Z (xy 



i 

e6£(l - 6?) 



(3.6) 



with 



4£ 



we identify as a fast- roll inflation (or rapid- roll) [64-66]. The first reason is that w c s calculated at 
this critical point can be made close to —1 especially for the phantom scalar field, and the second one 
is that the first coordinate of this point, using transformations (2.8) can be put in the following form 

<j) = -6£H<j> 

which for the conformal coupling £ = 1/6, reduces to condition for the rapid- roll inflation given by 
Kofman and Mukohyama in [65] . That is we identify this critical point as a generalisation to the non- 
minimally coupled case (both for the canonical and phantom scalar fields) with additional presence 
of the barotropic matter with the equation of state parameter w m . 
The linearization matrix for this critical point is the following 



\ 





dy' I 
dx I 2a 





U OX I 2a 

12£ 



(3.7) 



where 



dx' I 
d\ 1 2a 
dy' I 
dx I 2a 
dy' I 
OX I 2a 



-3(y* 2a ) 2 {l + w m + ^(l - 3£)A* 2 ^(Ay)/(AL), 

-el2£(l-6£)2/ 2 V(AU 

[(%*JMAL) + (e6(l - 6£)(1 + w m )z(X* 2a ) 



(V*2a) 2 )Ka )AKa) 



The eigenvalues of the linearization matrix are obviously l\ = 0, l 2 — 12£ and 1% 



-12£. Thus the 



fixed point is a non-hyperbolic and we cannot make any conclusions concerning its stability based 
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Figure 2. The phase portrait for the system (3.9) on the invariant submanifold in the vicinity of the critical 
point corresponding to the rapid-roll inflation. The bold line represents the center manifold for the problem. 
On the left diagram we present an unstable case for e(l — 3w m ) < and the right diagram for a stable case 
for e(l — 3w m ) > 0. The example is given for z(X) = £ and a = 1, w m = left for e = —1 and £ = 1/4, right 
for e = 1 and £ = 1/8. 

on linearization and the Hartman-Grobman theorem is not applicable [67, 68]. The answer to the 
question of stability or instability lies in the center manifold theory (see appendix A). 

We apply following procedure: first, we expand the right hand side of the dynamical system 
(2.17) into the Taylor series around the critical point (3.6) up to second order, and second, we make 
following change of dynamical variables 



= Po 



u 2a 
V - Via 

A - A; 



2« 



where the matrix P^ a is constructed from eigenvectors of the linearization matrix (3.7) calculated for 
corresponding eigenvalues and its inverse is 



P- 1 - 



1 

J_iV| i _ 

12$ dx I 2a 




1 dx' I 
12$ dX 1 2a 



1 ( 1 dx' I dy' I , dy' I \ 

24$ ^ 12$ d\ 1 2a dx 1 2a dX \ 2aJ 



1 



/ 



Then dynamical system in the vicinity of the critical point representing the fast-roll inflation is in the 
following form 

U 1 = -3j/2a( 2 C^2a z (^2a) + 1 + M m)«l' + A U VJ 2 + B U VW + C U UW, (3.8a) 

v> = -12£v + e±(l - 3w m )y* 2a u 2 - H {2^ a z(X* 2a ) + 1 + w m )y* 2a v 2 + e 12£(l - 60z(\L)uv 

+A v w 2 + B v vw + C v uw, (3.8b) 
w 1 = 12£w + A w w 2 + B w uw, (3.8c) 

where Ai , Bi and Ci are coefficients consisting of second derivatives of right-hand sides of dynamical 
system (2.17) calculated at the critical point under considerations. We can note that this dynamical 
system admits the invariant submanifold w = 0, and the dynamics can be well approximated on this 
submanifold. Then 

;, (3.9a) 



u = -3y2 a [2^\* 2a z(X* 2a ) + 1 + Wr, 
v' = -12£v + e^(l - 3w m )y* 2a u 2 ■ 



(2£\* 2a z(\* 2a ) + 1 + w m )y* 2a v 2 + el2f(l - 60^(A|Ju((?.9b) 
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on the invariant submanifold w = 0. 

From the center manifold theorem (appendix A) we have 

v = h(u) = - 3w m )y* a u 2 + ^^(1 - Su^MA^)^ 3 + 0(u 4 ) 

and inserting this approximation into (3.9a) we receive that the vector field restricted to the center 
manifold is given by 

i = -e^(l-3w m ) V 3 + 0( V 4 ), 

which indicates that for e(l — 3w m ) < it is an unstable and for e(l — 3w m ) > it is a stable critical 
point on the invariant submanifold w = (see figure 2 example for z(X) = — ). This equation can be 
simply integrated resulting in 

v{T f = 

/l ' (jp*)*e(l - 3w m )r + 1 

Above equation describes behaviour of the system on the center manifold which constitutes the in- 
variant submanifold. This solution can be used in construction of exact solution of the system (3.9) 
in the vicinity of the critical point representing the fast-roll inflation epoch. 

Using the solution from the center manifold theorem and keeping linear term in w only 

u(t)v(t) oc u(t) 3 rj 0, v(t)w(t) cx u(t) 2 w(t) « 0, v(t) 2 oc m(t) 4 rj 0, w(t) 2 rj 0, u(t)w(t) ps 

from (2.18) and (2.20) we get the parametric equations for the evolution of the scale factor and 
Hubble's function 

fln(^) = -e(l-6eMA^)z'(ALK"(cxp(12er)-l), 

* ; 2 . ., N V (3.10) 



(ln| 


V a 2a ) 






[ln( 







4? 

which can be easy combine resulting in 



Mw ini ( cxp(12£r) - 1 



H V A , / a 



ln [h™) e4Z(l - 60z(\* 2a )z>(\* a ) ln (3 - 11} 
where 

A = 7&~, n ^{ [ -el44f(l - 6^(1 + w m )(l + 3w m )z(\* 2a ) + 96^(2 - 9£)A| + 

+288e 2 (A^) 2 z(A^)]z'(A^J - 288£ 2 z(AL)}. 

One can conclude that needs \A\ <§; 1 in order to achieve iJ 2 « const, during the evolution. The 
linearised solution in w direction is valid up to the Lyapunov time r en d = using this we obtain 
maximal values of the scale factor and the Hubble's function respectively: 



ln 

V a 2a J 

-('fir) =~%A. 



:(l-6e>(A5> , (A$ )w fai (e-l) > 



(3.12) 



In figure 3 we present the evolution of Hubble's function and Hubble's horizon in the vicinity of this 
critical point. 

3.2.2 Slow-roll inflation 

The critical point located at 

* _ n ( * \2 2g(l - 3w m ) 2 _ 1 ,„ 

(1-60 { 2 &2b Z ( X 2b) + 0- + W m)) £6 ^ i 1 - 6 
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Figure 3. Evolution of InH 2 (left panel) and Rh (right panel) as a function of natural logarithm of the scale 
factor In a for a sample trajectory with e = — 1, £ = 4, z'(A*) = ^ = 10 in the vicinity of the critical point 
corresponding to the fast-roll inflation. The solid black line represents the linearised solution (3.10) and the 
dotted line represents the numerical solution of the system (2.17). 



with 



6£ 



wc identify as representing the phase of a slow-roll inflation due to x cx 
vicinity of this point corresponds to the slow-roll condition ^RiO. 
The linearization matrix is in the form 



A 



2 b 



dy' 

t).r 



dx' 




dx' 


2b dy 


21 1 


dX 


dy' 




dy' 


2b dy 


2b 


d\ 











so the dynamics in the 



(3.14) 



where nonzero elements are 

-$£_(! -3t« m ), 



dx' I 
dx 1 26 
dx I 
dy 1 26 

dx' I 
d\ 1 26 
dy' I 
dx 1 26 
dy' I 
dy 1 26 
dy' I 
9 A 1 26 



l-6f ' 



-18£(2£A* 6 2(A| 6 ) + (1 + ™ m ))j4*(A* 6 ), 
I^fe) 2 ( - e£ + (1 - 6£)z'(A^ b )(6£ 2 A* b z(A; b ) + (1 + «,„,))), 
e6f (1 - 3w m )y* b z(A^), 
—^(1-3^), 



1-65' 



-3&&((l&)MA») + ^'(^ b )(fe) 2 A^ + s6z(A* b )(l + w m - 8£)) 
and additionally we have the following relation 

6£ 



dx 



26 dy 



21, 



1 - 6£ 



(i - 3w m y 



The characteristic equation for the linearization matrix gives us vanishing eigenvalues l\ = I2 = I3 = 0, 
so the critical point is degenerated. In this case we cannot use standard procedures, following the 
Hartman-Grobman theorem [67, 68] of determining qualitative behaviour of the investigated system 
in the vicinity of this critical point. Instead we can notice that the linearization matrix A 2 b calculated 
at this critical point is nilpotent of order 3, i.e. (A2&) 3 = 0. Then solution of the linearised problem 
can be presented in the following form 



x(r) 



1 + A 2b T 
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Figure 4. Evolution of In H 2 (left panel) and Rh (right panel) as a function of the natural logarithm of the 
scale factor lna for a sample trajectory with e = +1, £ = |, z'{\\) = — = 100 in the vicinity of the critical 
point corresponding to the slow-roll inflation. The solid black line represents the linearised solution (3.15) 
and the dotted line represents the numerical solution of the system (2.17). 



Finally solutions in the vicinity of this degenerated critical point up to linear terms are 



where 



2b 



' 2b 



b 2b/ 



dx' 
dy 



( dx' 


dx' 




dx 


dy' 




I dx 


2b 9X 


H 

26 


- dy 


26 dX 


26/ 



I (^2b ~ AjJt 2 , 



y(r) = V2t + [ & Mfb - *2b) + ir M ~ *&) + °A J x * - *5*) 



dy' 



dy' 



2b 



2b 



26 



( dy' 


dy' 


H 


_ dy' 


dx' 




\ dy 


26 9X 


26 


dx 


26 dX 


26 / 



A(r) = Xf b \ 



dx' 


26 + 


dx' 


dy' 




26 d\ 


dy 


26 d\ 


26 


dy' 


26 + 




dx 




26 dX 


dx 


26 d\ 


26 



(3.15) 



36£(1 + w m ){2 + ^X* b z(X* 2b )) (y% b ) 2 z'(X* 2b ), 
18£(1 + w m )(\* b + £ 4(1 - 6£MA^))(2/* b )V(A* h ). 

These linearised solutions are valid up to a maximal value of the time parameter t = r max which can 
be used to calculate the scale factor growth during the slow roll inflation 



„,cnd 



In- 



The direct application of the linearised solutions (3.15) to (2.18) and (2.20) gives us the approx- 
imated evolution of Hubble's function in the vicinity of the critical point representing the slow-roll 
inflation (figure 4). 

3.3 Radiation domination epoch generated by non-minimal coupling 

Following critical point located at 

1 



x* 3a :g(x)=0, y* 3a = 0, \% a : z{\) 1 



e6£(l - 60 



(3.16) 



where g{x) = e(l — 4£ — w m )x 2 + e4£(l — 3w m )z(Xl a )x + 1 2 |^ (1 — 3w m ) and at this point value of 
the effective equation of the state parameter is 



WcK 
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represents the radiation dominated universe. With solutions to g(x) = equation in th form 



which is real only if the expression in square root is positive and for the barotropic matter with 
w m < | it is possible only if e£ < 0. We are interested only in evolution with £ > because of the 
discussion of the critical point representing the finite scale factor singularity, and this is the reason we 
identify this critical point as a representing radiation dominated epoch only for the phantom scalar 
field. 

The linearization matrix calculated at this point is in the following form 

dx I 

dx 1 3a " OA 1 3a 

A 3a =\ , (3.17) 




where 



dx' I 

d\ ha ■ 



^| 3a = £ 12£(1 -6Z)z(\* 3a )x* 3a , 

f£| 3a = e6£z'(Ay [2(1 - &0z(\* 3a )x* 3a + (1 - 3w m )(x* 3a + 6&(\* 3a ))(x* 3a + z(\* 3a )) 
f^| 3a = -el2£(l-6l;)z(\* 3a )x* 3a . 

Eigenvalues of the linearization matrix are l\ = — el2£(l — 6£)z(A 3a )x 3a , 1% = 0, l 3 = £l2£(l — 
6£)z(\ 3a )x 3a . This indicates that the critical point is non-hyperbolic one and the standard lineariza- 
tion procedure will be inefficient and we need to proceed with the center manifold theorem (see 
appendix A) and the procedure described during the discussion of the critical point representing 
fast-roll inflation. We make following change of dynamical variables 

x ~ x 3a 
y ~ Via 
^ - . 

where matrix P 3a is constructed from eigenvectors of the linearization matrix (3.17) and its inverse is 

P3- 1 = I 1 I , and X 

Then dynamical system can be presented in the following form 

v! = -£l2£(l - 6£)z{\* 3a )x* 3a u + A u u 2 + B u uw, 

v 1 = A v uv + B v vw, (3.18) 







dx' 






dX 


3a 






dX' 1 


<9x 


3a 


9 A ba 



el2£(l - 6^)z(Xl a )xl a w + A w u 2 + B w v 2 + C w w' 2 + D w uw, 



where Ai, Bi, Ci and Di are coefficients consisting of second derivatives of the right hand sides of 
dynamical system (2.17) calculated at the critical point (3.16). 

One can note that above dynamical system admits two invariant submanifolds namely v = and 
it = 0. 

On the first invariant submanifold the system can be simply reduced to 

u> = -el2£(l - 6^)z(X* 3a )x* 3a u . s 

w> = emi-60z(X* 3a )x* 3a w, [6 y > 

resulting in the solution representing a saddle type critical point in the form 

u(r) = u ini cxp ( - el2£(l - ^)z(X* 3a )x* 3a r) , 

w(t) = w™ exp ( e 12£(l - 60z(KaMaT) ■ 1 ' 
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Figure 5. The phase portrait for the system (3.21) on the invariant submanifold u — in the vicinity of the 
critical point corresponding to the radiation dominated universe with w e g = The bold parabola shaped 
line represents the center submanifold for the problem. On the left diagram we present an unstable case and 
on the right diagram for a stable case. The example is given for z(X) = ^, s = — 1, w m = 0, £ = 1/4 and 
a = 1 (left) and a = -4 (right). 



On the other hand we can also restrict our system to the invariant submanifold defined by u = 

then 

v' = B v vw, 

w + B w v 2 + C w w 

3a 

and from the center manifold theorem (see appendix A) we have 



w' = ¥- 

OX 



i^n.J (3.21) 



« = Hv) = -^r-v 2 + -J^- 3 {2B V - C w y + 0(v 5 ) 

dx l 3a v dx 

and inserting this approximation into first equation of the system (3.21) we receive that the vector 
field restricted to the center manifold is given by 



/ B V B W 3 

dx \ 3a 



where 



dx l 3a 



•2(l-6^(Ay 2 



and this indicates that for fj/^" 1 < it is an unstable and for f^y 3 " > it is a stable critical point 

on the invariant submanifold u = (see figure 5 for an example for z(A) = — ). 

Now we are ready to present the evolution of Hubblc's function in the vicinity of this critical 
point. First, we use approximated solutions (3.20) on the invariant submanifold v = 0. We have 

111 fe) = ( 6XP ( " £l2 ^ 1 ~ ^) z ^a)^ a r) - l) , 

In^) 2 = Y-er ^+ 6 ffi"^ ini (cx P (el2C(l - 6^z(XM a r) l)- (3.22) 
- (4 + ifc) ( ex P ( - ^2^(1 - 60.(AL)^ a r) - l) 

On the other hand, using the solution from a center manifold and keeping only linear terms in u 
w(t) 2 oc v(t) a ps 0, u(t)w(t) oc u(t)v(t) 2 ps 0, u(r) 2 ps 
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Figure 6. Evolution of InH 2 (left panel) and Rh (right panel) for a sample trajectory with e = — 1, £ = 10, 
z '(^3a) = a = To m the vicinity of the critical point corresponding to the radiation dominated epoch for 
the phantom scalar field. The solid black line corresponds to the linearised solution (3.23), the dashed line 
corresponds to exact radiation dominated expansion of the universe InH 2 cx —4 In a and the dotted line 
corresponds to the numerical solution of the system (2.17). 



wc receive 



ln ( J* ) = ^u 1111 ( exp ( - el2f (1 - *Lt) - 1 



2 (3 23) 

(*) = ~( 4 + ^ z gi } ) ^u^(e X p(-emi-60z(XM a r) - l) 



In 

which can be easy to combine as 



One can notice that this expression resembles behaviour of the Hubble's function during the pure 
radiation domination epoch, but with contribution coming from non-minimal coupling. The linearised 
solutions are valid up to the Lyapunov time r en d = _ e i2g(i-6g)z(A* )x* > ^' t nen inserting this in to 
the latter equations we receive maximal values of the scale factor and Hubble's function valid in the 
center manifold approximation 




(.-.). 

4? ^AfTJ 



^Hu ini [e- 1 



(3.25) 



In figure 6 we present evolution of InH 2 and Rh as a functions of in a in the vicinity of the 
critical point representing radiation domination epoch for the phantom scalar field. 

There is another critical point which represents the radiation dominated universe located at 

a& = 0, 24 = 0, \* 3b : z(X) 2 = -L (3.26) 

with the effective equation of the state parameter 

v* = \ 

We need to stress that this critical point exists only if w m ^ |. Linearised solutions in the vicinity of 
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this critical point are 

+ —L-( X ™ - (1 - 3 Wm )z'(\* 3b )(\™ - X* 3b )) exp (Z 3 r), (3.27a) 
lft»(r) = l^exp&r), (3.27b) 
A 3fc (r) = + _J_-JL^ (4f + z'(A5 b )(A™ - A3 b )) exp (Z lT ) - 

?7Ft(^ - (1 - 3^ m )z'(A^ b )(A™ - A 36 )) cxp(/ 3 r). (3.27c) 

where 

h = 6f(l - 3to TO ), Z 2 = 12£, Z 3 = -6£ 

are the eigenvalues of the linearization matrix calculated at this critical point. Simple inspection of 
this eigenvalues gives us further constraint on the value of the barotropic matter equation of state 
parameter w m , namely l\ should be positive resulting in w m < | to assure that in (x,X) plane the 
dynamics in the vicinity of this critical point would correspond to a saddle type critical point. This will 
guarantee that the evolution proceeds towards the next critical point representing matter dominated 
universe. 

Using linearised solutions (3.27) we are able to express (2.18) and (2.20) as a parametric functions 
of time r 



' hl (#) = 6 ^ ~ g2 (i-s^jgU^) *(*&) (4? + z '(Kb) (A£ - K„) ) ( exp (hr) - 1) - 

(A3 6 )(4fc " (1 " 3 Wm )^'(A 3b )(A^ - A* 6 )) (exp(Z 3 r) - l), 

24^ - e 2^1 - (1 _ 3 :;Kfi„,) ) z (^)(^ + a'C^K^ - K b )) (exp far) - l)- 



, 2 

HK 1 



+ e8 2 i ^7^ A 5 6 )(4" i - (1 - 3«; m )z'(A5 b )(A 3 ' b i - AJ 6 )) (cxp(Z 3 r) - l) 

(3.28) 

The linearised solutions (3.27) are valid up to the Lyapunov characteristic time r cn d = j~ = ilj ancl 
inserting it in the latter equations we can obtain maximal values of the scale factor and the Hubblc's 
function valid in the linear approximation. 

The zero-order approximation (.t™ = 0, A™ = A5 b , yfj ^ but (y™) 2 0) is 

(3.29) 




and it can be combine to 

' a 



which is the exact behaviour of Hubblc's function during the radiation domination era. This approx- 
imation is also valid up to T cn , 
the scale factor grows at least 



imation is also valid up to r cn d = so one can calculate that during radiation domination epoch 



end _ mi f~ 

a 3b - a 3b Ve- 
in figure 7 we present evolution of In H 2 and Rh as a function of In a in the vicinity of the critical 
point representing the radiation domination era. 

3.4 Matter domination 

The next critical point is located at 

x* 4 = 0, y% = Q, \l:z{X) = Q (3.30) 
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Figure 7. Evolution of InH 2 (left panel) and Rh (right panel) for a sample trajectory with e = +1, £ = y^, 
z '(^sb) = a ~ 100 m the vicinity of the critical point representing the radiation dominated epoch for the 
canonical scalar field. The solid black line corresponds to the linearised solution (3.28), the dashed line 
corresponds to exact radiation dominated expansion In if oc —4 In a and the dotted line corresponds to the 
numerical solution of the system (2.17). 



and w a ff given by (2.12) calculated at this point is 

W c g = W rn . 

We identify this critical point as representing the universe which dynamics is dominated by the 
barotropic matter included in the model with the equation of state parameter w m . 
The linearised solutions in the vicinity of this critical point are in the form 

X4 ( T ) = A K ni - z'(AI)k(Af - Al)) exp (I lT )- 

- A( x 4 ni " z'(A2)Zi(Af - AI)) exp (1 3 t), 
y 4 (r) = </f exp(Z 2 T), (3.31) 
A4 ( r ) = A 4 + H - z '( X *M^f ~ K)) cx P (Jit)- 

2'(AJ)^l-( 3 J V / 

xr - z'(AI)?i(Ar - AI)J exp (/ 3 t). 



2'(A2)Ui-i 3 



where 



h = -\ ki - «>m) + y(i - ™ m y - ye(i - 3w m ) ) . 

3 

h = + 

k = ~ ^(1 - » m ) - \j (1 - «->m) 2 - y^ 1 _ 3u, m) 

We need to note that this critical point can became degenerate for two specific values of w rn , namely, 
for w m = — 1 the second eigenvalue vanishes and for w m = | the third eigenvalue vanish, for any value 
of the coupling constant £, which makes the system in the vicinity of this critical point degenerated. 
From linearised solution (3.31) we have 

x 4 (r) 2 «0, 2/ 4 (t) 2 «0, z(A 4 (t)) 2 «0 

and from (2.18) and (2.20) parametric equations for evolution of the scale factor and Hubble's function 
are 

t, 

(3.32) 

-3(1 + W m )T. 



fin) 








[H 
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Figure 8. Evolution of InH 2 (left panel) and Rh (right panel) for a sample trajectory with e = —1, £ = i, 
z'(Ai) = — = 100 in the vicinity of the critical point corresponding to the barotropic matter dominated 
universe. The solid black line represents the linear approximation (3.32) and the dotted line represents the 
numerical solution of the system (2.17). 



Combining these two expressions we get the Hubble's function as function of the scale factor during 
the barotropic matter domination epoch 

-3(1+™,, 

H 2 = (Hf) 2 ' 

The linearised solutions (3.31) are valid up to the Lyapunov time r en d = 77 

..(fl) - • 



.mi 

M 



3(l+w m ) ' 



hi 



(Hf 
\ HI 



3(l + i 

= -2. 



(3.33) 



One can notice that for dust matter w m = during the matter domination epoch the scale factor at 

a ond 2 

least grows -^r = d w 1.948 times. 

In figure 8 we present evolution of the scale factor and Hubble's function in the vicinity of the 
critical point representing the barotropic matter domination epoch. 

3.5 The present accelerated expansion epoch 

Finally we proceed to the last critical points located at 



4 = 0, (y* 5 ) 2 = l-e6£z(\* 5 ) 2 , A^(A) 2 +4z(A)- — =0 



e6£ 



(3.34) 



with 



W e f{ = -1. 



There can be more than one such critical points because of the third equation in (3.34) which can 
have more than one solution. In what follows we will show that at least one of them represents a 
stable critical point. 

In this case the characteristic equation for eigenvalues of the linearization matrix calculated at 
this critical point is in the form 

I 3 +pl 2 + ql + r = 

where 



q = (l-e6£(l-6£MA^) 2 ) 



r = 3(l + m m )(l-e6£(l-6£MA*) 2 ) 2 



12£(1 + e6^(A^) 2 ) + 9(1 + m m )(l - e6£(l - 6£)*(A^) 2 ) 



1 fa.;) 4 

2 z'(\;j 



12£(1 + e6^(A* 



In the most general case without assuming any specific form of the potential function we are 
unable to solve this equation. In spite of this we are able to formulate general conditions for stability 
of this critical point. This requires that the real parts of the eigenvalues must be negative. From the 
Routh-Hurwitz test [68] we have that the following conditions should be fulfilled to assure stability 
of this critical point 



p > 



and 



r > 



and 



q-->0 
P 



Simple inspection of these conditions gives us that, for any matter with w rn > —2, p is always positive 
because of (l — e6£(l — 6£)z(Xl) 2 ) > due to time transformation and that if r is a positive quantity 



it follows that q 
critical point 



is positive too. We conclude that if the following condition is fulfilled at the 



Re[h,2, 3 ] < 



— £ 



12£(1 + e66?(A£) 2 ) >0 



■2z'(X* 5 ) 

it represents a stable critical point with the negative real parts of the eigenvalues. 
In order to simplify this condition let us introduce the following function 



(3.35) 



h(X) = \z(X) 2 + 4z(A) 



A 

£6£' 



(3.36) 



where location of the critical point is the solution to the equation h(X) = 0, and obviously h(X£) = 0. 

Let us assume that A5 ^ it follows from (3.36) that also z(X§) ^ but h(X$) = 0. Differentiation 
of Eq. (3.36) gives 

h'(Xt) = z{Xlf - ^ + 2z'{Xt){Xlz{Xl) + 2) 
which after little algebra can be transformed in to the following form 



£6£ 



l+e6£z(X* 5 y 



and finally we arrive to the reformulated stability condition (3.35) in the form 



Re[h,2,z] < 



> o 



(3.37) 



We have reduced analysis of stability of the critical point representing accelerated expansion to the 
simple analysis of the sign of the quantity given by relation (3.37). If we assume that the function 
z(X) is a monotonic one, i.e. it is a growing or decreasing function in the interesting region of the 
phase space then it follows that if we have at least two critical points given by (3.34) one of them is 
definitely a stable critical point. 

To present the evolution of Hubblc's function in the vicinity of this critical point, as an example, 
we choose the simple form of 2 (A) = — function, and values of the parameters £ and a in range 
for which there exists only one critical point corresponding the present accelerated expansion of the 
universe [54]. The linearised solutions in the vicinity of the critical point located at x§ = 0, (1/5) 2 = 1, 
A? = are 



(3 + VSi) [4 ni + 5s(3 - V^)Xr\ exp (/it) 

-(3 - 7A^) [4- + 5^(3 + v^Af ] exp (l 3 r) 
2/5 (r) = yt + (yr-yt)exp(l 2 r), 



(3.38) 



Ab(t) 



ini 1 
C 5 + 2^ 



(3-V&E)Xf expftr) 



x 5 



^(3 + v^)Af exp(/ 3 r) 



- 19 - 



Figure 9. Evolution of InH 2 (left panel) and Rh (right panel) as a function of the natural logarithm of 
the scale factor lna for a sample trajectory with e = — 1, £ = 1, z'(\l) = ~ = 20, j/5 111 — yt = — j^j in the 
vicinity of the critical point representing the present accelerated expansion of the universe. The solid black 
line represents the linear approximation (3.39) and the dotted line corresponds to numerical solution of the 
dynamical system (2.17). 

where l h3 = -^(3± ^9 + s2a - 48£) and l 2 = -3(l + w m ) are eigenvalues of the linearization matrix 
and A 5 = 9 + e2a - 48£. 

Then keeping only linear terms in initial conditions 

x 5 (r) 2 «0, z(A 5 (t)) 2 «0, A 5 (t)z(A 5 (t)) w 0, 

y 5 (T) 2 « {ytf + 2y* 5 (yf - y%) exp (l 2 r), 

from (2.18) and (2.20) we receive the parametric equations of evolution of the scale factor and Hubble's 
function 



flnj 


V a T ) 






[in 


:+) 





/ n (3-39) 

In = 2y* 5 (yf - y*) (l - exp ( - 3(1 + w m )r)J 

Combining these two expressions we get Hubble's function as a function of the scale factor in the 
vicinity of the critical point corresponding to the present accelerated expansion of the universe 



H = (HD 2 exp i 2y* b (y™ - y* 5 ) ( 1 - ^ 



-3(l+w m ) 



One can notice that taking the following limit 

H 2 Rn = lim H 2 = Urn H 2 = (H™) 2 exp {2y* 5 (y™ - y* 5 ) } « (H™) 2 (l + 2y* (yj* - y* 5 )) 

we get the asymptotic de Sitter expansion. 

In figure 9 we present the evolution of In H 2 and Rr as a function of In a in the vicinity of this 
critical point. 

4 Summary and Conclusions 

Modern cosmology becomes very similar to the particle physics. Both theories have parameters and 
characteristic energetic cut offs. They are the effective description of deeper physics which is currently 
unknown. The values of these parameters should be obtained form more fundamental theories or from 
observations. In cosmology (ACDM model is called Standard Cosmological Model) the role of such 
a parameter plays the cosmological constant. Our proposition is to extend this paradigm in which 
matter content is described in terms of barotropic perfect fluid by introduction additional scalar field 
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Figure 10. The phase space portrait for the model with the cosmological constant and the canonical scalar 
field (e = +1) with £ = 1/8 and the dust matter w m = 0. The critical points are: S - the finite scale factor 
singularity, RI - the rapid-roll inflation, SI - the slow-roll inflation, R - the radiation dominated era, M - 
the barotropic matter dominated era and Q - the quintessence era. Note that the critical points representing 
the finite scale factor singularity, the rapid-roll inflation and the slow-roll inflation have the same value of 
coordinate z. 

non-minimally coupled to gravity. As a result we discover new evolutional path which open new 
perspectives of description of cosmological evolution in unified way. In this scheme the inflation era 
appears in natural way and it is not put into the ACDM scenario by hand. 

In this paper we have shown that the all important epochs in the evolution of the universe can 
be represented by the critical points of the dynamical system arising from the non-minimally coupled 
scalar field cosmology in spite of not assuming a form of the potential function. We have shown that 
for the positive coupling constant there exists a past finite scale factor singularity for both types of the 
scalar fields. Additionally all the intermediate states are transient one, i.e. they are represented by 
an unstable critical points in the phase space and last for an finite amount of time. The existence of 
the radiation dominated era is purely the result of the evolution of the non-minimally coupled scalar 
field. 

For the canonical scalar field e = +1 and 0<£<l/6we can construct the unique evolutional 
path represented by the trajectory in the phase space which travels in the vicinity of the following 
critical points (figure 10) 

and for the phantom scalar field e = — 1 and £ > 1/6 (figure 11) 

ln>2ai->-3ai->4i->5 
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Figure 11. The phase space portrait for the model with the cosmological constant and the phantom scalar 
field (e = — f) with £ = 1/4 and the dust matter w m = 0. The critical points are: S - the finite scale 
factor singularity, RI - the rapid-roll inflation, R - the radiation dominated era, M - the barotropic matter 
dominated era and Q - the quintessence era. In the case of the phantom scalar field the critical point 
representing slow-roll inflation is not present. The critical pints denoted as 5", RI and R have the same value 
of coordinate z. 



Within one framework of non-minimally coupled scalar field cosmology we were able to unify 
all the major epochs in the history of the universe (see figure 12 for twister type behaviour where 
trajectories interpolate between radiation era, matter domination era and quintessence epoch). 

From the analysis presented in this paper one can draw the general conclusion that if the non- 
minimal coupling constant is present and is different from the conformal coupling £ ^ 1/6 then new 
evolutional types emerge forming the structure of the phase space nontrivial and richer. Moreover the 
coupling constant gives us the effect of continuation (glues the evolution) between different cosmo- 
logical epochs which is very attractive in cosmology, serving as a potential explanation of the global 
properties of the universe. 

A The Center Manifold Theorem for three-dimensional dynamical sys- 
tems 

For the sake of completeness we present here the theorem concerning behaviour of nonlinear dynamical 
system in the vicinity of degenerated critical point. Expanded discussion can be found, for example, 
in books by Perko [67] or Wiggins [68]. 

Suppose we consider the following 3-dimensional nonlinear dynamical system 

x = /(x), xel 3 (A.l) 
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Figure 12. The phase space portrait representing twister type behaviour. Trajectories in this type solution 
interpolate between three major epochs in the history of universe: R - the radiation dominated universe 
with U! e ff = ^, the matter domination epoch (an unstable focus type critical point) and Q - the quintessence 
domination epoch with w e ff = — 1. This type of evolution does not depend on the form of assumed function 
2(A) (i.e. the form of the scalar field potential) and is generic for the canonical scalar field cosmologies (e = +1) 
with £ > and the barotropic matter with equation of state parameter — 1 < w m < 1/3. 

We are interested in the nature of solution to this dynamical system near fixed point x for which 
/(x) = 0. 

First, we transform the fixed point x = x of (A.l) to the origin using the transformation y = x— x. 
Then the system (A.l) becomes 

y = /(x + y), yeM 3 (A.2) 
then Taylor expansion of /(x + y) about x = x gives 

y = ,4y + i?(y), y g R 3 (A.3) 

where A = Df(5t) ia a linearization matrix calculated at the fixed point, R(y) = 0(\y\ n ) and we have 
used /(x) = 0. 

From now on we will assume that the linearization matrix has purely real eigenvalues and one is 
zero l\ = 0, one positive I2 > and one negative I3 < 0. From elementary linear algebra we can find 
a linear transformation P which transforms the linear part of equation (A.3) into a diagonal form 

(A.4) 
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with a linear transformation of variables 



and the matrix P is constructed from the corresponding eigenvectors of the linearization matrix 
A. Using this same linear transformation to transform the coordinates of the nonlinear part of the 
system (A. 3) gives the following 

u = F(u,v,w), 

v = l 2 v + G[u, v, w) , (A. 5) 

ill = I3W + H(u,v,w\. 

where F(u,v,w), G(u,v,w) and H(u,v,w) are polynomial in the coordinates. The fixed point 
(it, v, w) = (0, 0, 0) is unstable due to the existence of a 1-dimensional unstable manifold associated 
with the negative eigenvalue I3 < 0. 

Definition A.l (Center Manifold) An invariant manifold will be called a center manifold for (A. 5) 
if it can be locally represented by 

W c (0) = \{v,,v,w) e R 3 \v = h±(u),w = h 2 (u),\u\ < 5,^(0) = 0,/i-(0) =0,i = 1,2} (A.6) 

for S sufficiently small. 

Theorem A.l (Existence) There exists a C r center manifold for (A. 5). The dynamics of (A. 5) 
restricted to the center manifold is, forr/ sufficiently small, given by the following 1-dimensional vector 
field 

f) = F(r,,h 1 (r,),h 2 (r,)). (A.7) 

From the fact that the center manifold is invariant under the dynamics generated by (A. 5) we 
obtain 

u = F(u, hi(u), h 2 (u)), 

v = h[(u)u = hhi(u) + G(u,hi(u),h 2 {u)) , (A. 8) 

w = h' 2 (u)u — hh 2 (u) + ij(u, /ii(u), /i 2 (u)) , 

which yields the following quasilincar differential equation for h\{u) and h 2 (u) 

Af(hi(u)) = /ii(u)F(u, In(u), h 2 (u)) - l 2 h!(u) - G(u, hi(u),h 2 (u)) = 0, , 
M(h 2 (u)) = h' 2 (u)F(u, hi(u),ha(v)) - hh 2 {u) - H(u, hi{u), h 2 {u)) = 0, 1 y) 

and the following theorem (see Theorem 18.1.4 in Wiggins [68, p. 248]) justify solving (A. 9) approx- 
imately via power series expansion: 

Theorem A. 2 (Approximation) Let <f>: K — > K be a C 1 mapping with 0(0) = (j>'(0) = such that 
Af(4>(u)) = 0(\u\ q ) as u for some q > 1. Then 

\h(u) - cj)(u)\ =Q(\u\ q ) as u^0. 

This theorem allows us to compute the center manifold to any desired degree of accuracy by solving 
(A. 9) to the same degree of accuracy. 
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